Discovery of a novel bat lyssavirus in a Long-fingered bat (Myotis capaccinii) from Slovenia

Lyssaviruses are the causative agents of rabies, a zoonotic, fatal disease that is thought to be ancestral to bats. In the last decade, the detection of bat associated lyssaviruses is increasing also in Europe. Within a retrospective bat associated lyssavirus surveillance study a total of 225 dead bats of 21 bat species were collected in Slovenia between 2012 and 2019 and tested by specific real-time RT-PCR method. The first lyssavirus positive sample in bats in Slovenia was detected using the real-time RT-PCR, the fluorescent antibody test, and next generation sequencing, while the rabies tissue culture inoculation test was unsuccessful due to sample degradation and storage conditions. The nearly complete genome of Divača bat lyssavirus from Slovenia consists of 11,871 nucleotides and reflects the characteristic gene organization known for lyssaviruses, encoding the five viral proteins. Phylogenetic analysis of Divača bat lyssavirus revealed that it belongs to phylogroup I lyssaviruses and is most closely related to Kotalahti bat lyssavirus (KBLV) with 87.20% nucleotide and 99.22% amino acid identity. Together with KBLV, Khujand virus, European bat lyssavirus 2, Bakeloh bat lyssavirus, and Aravan virus, Divača bat lyssavirus was detected in the genus Myotis suggesting its key role in the transmission and maintenance of certain lyssaviruses.


Introduction
Bats belong to the mammalian order Chiroptera, with more than 1460 species in 21 bat families [1]. They are widely distributed on all continents except Antarctica and play an important role in seed dispersal, pollination, and in organic fertilization through their guano [2]. Bats are a natural reservoir for viruses that are among the deadliest viruses transmitted from wildlife to humans, including rabies, and SARS coronavirus [3]. As reservoir hosts for emerging viral zoonotic diseases, bats also carry viruses from at least 28 different virus families, most of which are likely host-specific and have limited zoonotic potential [4]. The number of viral zoonoses in bats is comparable to the number of viral zoonoses in other mammal or bird orders since variation in the number of viral zoonoses among animal groups arises as a consequence of their species richness [5].
In this manuscript, we describe the discovery and genetic characterisation of a previously unknown bat-associated lyssavirus, Divača bat lyssavirus, in an insectivorous bat species from Slovenia as a result of national retrospective surveillance programme. Based on whole genome sequencing and phylogenetic analysis, Divača bat lyssavirus, together with KBLV, may represent a putative new lyssavirus species in Europe.

Sampling
Between 2012 and 2019, a total of 225 dead bats of 21 bat species from 59 (out of 212) municipalities in Slovenia were collected by bat biologists or volunteers and stored frozen (S1 Table). The collected dead bats were submitted to the National Veterinary Institute, Institute of Microbiology and Parasitology, Virology unit, for lyssavirus diagnosis. Bat biologists provided information on bat species, location (municipality), and year of collection. Bat species were identified using morphological keyes, described in Dietz and Halversen [35]. Additional molecular characterization of the host species from lyssavirus positive bat sample PP-0868/ 2014, where NGS reads were mapped to cytochrome b and cytochrome c oxidase subunit I [36], confirmed the original assignment classification to Myotis capaccinii.
Brain samples were collected through the foramen occipitale magnum by pipette aspiration. After aspiration of the brain tissue, the cranial cavity was rinsed with RPMI-1640 medium (Thermo Fisher, USA). The collected brain tissue was homogenized in a total volume of 500 μl of RPMI-1640 medium (Thermo Fisher, USA) before long-term storage at < -60˚C.
Permit for capture, disturbance, and temporary taking from the wild and sampling of protected animals was issued by the Agency of the Republic of Slovenia for the Environment No. 35601-35/2010-6.

Diagnostic methods
For molecular detection of lyssaviruses, brain tissue homogenates were subjected to automated RNA extraction using the KingFisher Flex Purification System (Thermo Fisher, USA) and MagMax Core Nucleic Acid Purification Kit (Thermo Fisher, USA) according to the manufacturer's instructions. The extracted RNA was tested by real time RT-PCR [37,38] to detect lyssaviruses. Real-time RT-PCR with three specific primers LN34forward1, LN34forward2, and LN34reverse, and two probes LN34probe and LN34probeLago targeting highly conserved sequences within the nucleoprotein gene region [37] was performed.
When the lyssavirus-positive sample PP-0868/2014 was detected by real-time RT-PCR, the bat carcass was resampled by opening the cranial cavity to obtain the brain sample, which was subjected to a fluorescent antibody test (FAT) performed as previously described [39].
RTCIT was performed from lyssavirus positive sample PP-0868/2014 as previously described [39]. Three consecutive serial passages were performed to obtain RTCIT results.

Whole genome sequencing and phylogenetic analysis
The lyssavirus positive sample PP-0868/2014 was subjected to next generation sequencing (NGS) to determine the complete genome sequence of lyssavirus using the NGS service of Novogene (Cambridge, UK). From extracted RNA to final data, Novogene's service included sample preparation, quality control, library construction using the RIP-seq library preparation kit (Illumina, USA), library quality control, sequencing on NovaSeq 6000 (Illumina, USA), using the NovaSeq PE150 kit (Illumina, USA), and data quality control. NGS reads were used for de novo assembly with SPAdes 3.13.0 [40]. Diamond BLASTx [41] and MEGAN 6.11.7 [42] were used for the taxonomic assignment of the assembled contigs. The nucleotide sequence of the assembled lyssavirus genome, named Divača bat lyssavirus (original name of sample PP-0868/2014), was deposited in the GenBank database with accession number OQ428158 and compared with other genomes of 17 lyssavirus species and 2 putative lyssavirus species and annotated using Geneious 20221.1 (Biomatters, New Zealand). The nucleotide and amino acid alignments were constructed using MAFFT [43]. Phylogenetic analysis was performed using IQ-TREE 1.6.12 [44], with ModelFinder [45] determining the best model according to BIC score and 1000 ultrafast bootstrap replicates [46,47] used to test tree reliability.

Virus detection by real-time RT-PCR, FAT, and RTCIT
One sample (original name PP-0868/2014) out of 225 tested samples from 21 different bat species (S1 Table) was found to be lyssavirus positive by the real-time RT-PCR, with a cycle threshold value of 21.82. The sample PP-0868/2014 was from species Long-fingered bat (Myotis capaccinii), which was found dead approximately 700m inside of Škocjan caves in Divača municipality collected in the year 2014. The bat carcass was found in the advance stage of decay and was stored at < -20˚C until brain collection in 2020. The lyssavirus positive brain sample PP-0868/2014 was subjected to further analysis, namely FAT and RTCIT. The FAT slide stained with FITC anti-rabies monoclonal globulin showed positive results with typical lyssavirus staining (S1 Fig). No viral growth was detected on RTCIT for PP-0868/2014 after three consecutive passages.

Whole genome sequencing and phylogenetic analysis
Viral RNA from lyssavirus positive sample PP-0868/2014 was subjected to NGS to determine the complete genome sequence, named Divača bat lyssavirus. NGS yielded 13,048,710 raw reads that were used for de novo assembly. The nearly complete lyssavirus genome sequence with a total length of 11,871 nucleotides (an average sequence depth of 135,333 nucleotides per nucleotide site) was generated, encoding five genes, N, P, M, G, and L in the highly conserved order typical of lyssaviruses. The lengths of the N, P, M, G, and L protein genes were 1356 nt, 894 nt, 609 nt, 1581 nt, and 6384 nt, respectively. Gene lengths and intergenic spacers of Divača bat lyssavirus were found to be similar to all lyssaviruses from phylogroup I.
Phylogenetic analysis of the concatenated N+P+M+G+L coding sequences of the assembled Divača bat lyssavirus, 17 lyssavirus species, and 2 putative lyssavirus species (Fig 1) showed that Divača bat lyssavirus clustered with phylogroup I lyssaviruses (Fig 1) and is most closely related to KBLV with 87.20% nucleotide and 99.22% amino acid identity (Table 1). Similar nucleotide and amino acid identities to KBLV were found for the N gene and also for the P, M, G, and L genes. The percent nucleotide and amino acid identities of the N, P, M, G, L, and concatenated N+P+M+G+L genes between Divača bat lyssavirus and other lyssavirus species from GenBank are shown in Table 1.

Discussion
Several bat lyssavirus surveillance systems have been implemented in Europe and results have shown that passive or retrospective surveillance is better than active surveillance for lyssavirus detection [27,50]. In this study, a retrospective survey was conducted between 2012 and 2019, collecting 225 carcasses from 21 of the 32 bat species living in Slovenia [51]. Lyssavirus was found in one of the two Myotis capaccinii carcasses. Based on a single detection of lyssavirus in a particular bat species, we could not determine whether Myotis capaccinii is the true reservoir host for the Divača bat lyssavirus. Similar results and conclusions were reported by Nokireki et al. [34] who described the putative new lyssavirus species KBLV in Myotis brandtii and also found only one KBLV positive sample. However, the authors of a later study on KBLV [52] suggested that, based on the data that closely related phylogroup I lyssaviruses have also been isolated from bat species of the genus Myotis, Myotis brandtii is probably a reservoir host and not infected by spillover infection from another bat species. This study is the first to report a lyssavirus detected in Myotis capaccinii, a bat species in which several viruses from the families Astroviridae (3), Coronaviridae (6), Herpesviridae (2), and Paramyxoviridae (1) have been described [53]. Infection of bats with lyssaviruses occurs worldwide, although different virus species occur in different regions, have co-evolved and are therefore linked to specific bat species [54]. In the Americas, only RABV is associated with bats, whereas in Europe, Africa, Asia, and Australia, the other lyssaviruses predominate without RABV being associated with bats [55]. In Europe, numerous bat species have already been identified as infected with lyssaviruses [56]. However, there is no report of a lyssavirus finding for the bat species Myotis capaccinii.
Using NGS, a nearly complete genome sequence of the Divača bat lyssavirus was generated and phylogenetic analysis was performed. According to the ICTV [9], there are several criteria for is that a novel lyssavirus species occupies a distinct ecological niche, as evidenced by host species, pathobiological properties, or geographic range. The Divača bat lyssavirus occupies a specific ecological niche because it originates from a different host than KBLV, which has also been detected in a different geographical location (Fig 2). Taking all four ICTV criteria together to taxonomically classify the novel Divača bat lyssavirus, even though it was found in a new host species almost two thousand kilometres away from the closely related KBLV, we suggest that the Divača bat lyssavirus and KBLV are both members of the same putative novel lyssavirus species. As already aforementioned, Divača bat lyssavirus is most closely related to KBLV. KBLV was detected in Myotis brandtii in Finland in 2017 [34]. In Slovenia, Myotis brandtii is very rare, however, it was found in the vicinity of Škocjan caves [57]. The contact between Myotis capaccinii, species with circum-Mediterranean distribution [48] almost exclusively roosting in caves, and Myotis brandtii, species with wider distribution in Europe [49] primarily roosting in tree holes and cracks, remains unclear as they usually do not share roosts. The geographical distance between Kotalahti in Finland and Divača in Slovenia is quite large, it measures 1886 km (Fig 2). Due to the vast geographical distance and fact that both host species Myotis brandtii and Myotis capaccinii are mostly regional migrants [58], we could speculate that lyssaviruses similar to KBLV and Divača bat lyssavirus could be present in other European countries between Slovenia and Finland, suggesting the need for surveillance programs and the need to increase public awareness.
Despite efforts to isolate viable virus from the sample PP-0868/2014, we were unsuccessful and could not evaluate the vaccine efficacy of conventional rabies vaccines. However, according to the assignment of Divača bat lyssavirus to phylogroup I and high sequence identity with KBLV we could speculate good cross protection by conventional rabies vaccines.
We believe that the complete genome sequence of Divača bat lyssavirus and its phylogenetic relationships to other lyssaviruses identified in our study may be useful for future studies of the lyssavirus genus. For further studies, the isolation from field samples on cell culture will be an important step toward the possibility of the additional analysis of this strain.
In conclusion, Divača bat lyssavirus is most closely related to KBLV, which has been proposed as a separate species within the genus Lyssavirus. Previously published studies on KBLV indicate that it belongs to phylogroup I of the lyssaviruses against which rabies vaccines provide protection. Continued passive surveillance would be required to obtain infectious lyssaviruses for evaluating the efficacy of vaccines relevant to public health. Because Divača bat lyssavirus has the greatest genetic similarity to KBLV, it could be considered to pose a low risk to public health if preventive measures such as vaccination of bat caretakers and prophylactic treatment after human-bat contact are followed.  [48,49]. Base layer of the map was made with Natural Earth (https://www.naturalearthdata.com/http//www.naturalearthdata.com/download/10m/cultural/ne_10m_admin_0_sovereignty.zip).